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ABSTRACT 

The coexistence of Planck and Fermi satellites in orbit has enabled the ex- 
ploration of the connection between the (sub-) millimeter and 7-ray emission in 
a large sample of blazars. We find that the 7-ray emission and the (sub-) mm 
luminosities are correlated over five orders of magnitude, L 7 oc L( sub ^ mm . How- 
ever, this correlation is not significant at some frequency bands when simultane- 
ous observations are considered. The most significant statistical correlations, on 
the other hand, arise when observations are quasi-simultaneous within 2 months. 
Moreover, we find that sources with an approximate spectral turnover in the mid- 
dle of the mm-wave regime are more likely to be strong 7-ray emitters. These 
results suggest a physical relation between the newly injected plasma components 
in the jet and the high levels of 7-ray emission. 

Subject headings: galaxies: active; BL Lacertae objects: general; quasars: general; 
galaxies: jets; gamma rays: general; submillimeter; radio continuum:galaxies 
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Introduction 



Within the three years of its operation, the Fermi satellite has confirmed that 



the extragalactic 7-ray sky is dominated by emission from blazars (The Fermi-LAT 



Collaboration (2LAC) 2011 Ackermann et al. 2011b). A blazar is an unusual type of active 
galactic nuclei (AGN) in which a relativistic jet points toward Earth, and due to relativistic 
effects the emission is amplified and variability time scales look shorter. Blazars come in 
two flavours, with prominent or weak broad-emission lines, and are called Flat Spectrum 
Radio Quasars (FSRQ) or BL Lac objects (BLLacs), respectively. 

However, during the Fermi era it has been shown that not all blazars are strong 7-ray 



emitters (e.g. Lister et al. 2011 Giommi et al. 2011), and the likelihood of detection is 



probably related to faster and brighter jets (e.g. Kovalev et al. 2009 Nieppola et al. 2011), 



jets pointing closer to our line of sight (Lister et al. 2009), larger apparent jet opening 



angles (Pushkarev et al. 2009 Ojha et al. 2010), higher variability (e.g. Richards et al 



2011), or the presence of multiple inverse- Compton components (Abdo et al. 2010a Tiirler 



& Bjornsson 2011). In addition, it has been shown that the highest levels of 7-ray emission 



are closely related to ejections of superluminal components (e.g. Agudo et al. 2011) and 



ongoing high-frequency radio flares ( Leon-Tavares et al. 2011). Whether a blazar needs to 



fulfill all (or a specific combination) of the above conditions in order to radiate in 7-rays is 
one of the most important questions that still needs to be answered about the nature of 
7-ray emission in blazars. 

The complete diagnostic of what makes a blazar bright at 7-ray wavelengths requires 
measurements of the so far poorly explored (sub-)mm spectral bands. Because mm and 
sub-mm radiation in blazars merely samples the optically thin regime of the synchrotron 
spectrum, thermal emission (from accretion processes) and radiation from lobes is negligible 



in these bands (Giommi et al. 2009). Therefore, (sub-)mm measurements can efficiently 
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probe the inner regions of blazar jets and shed light on the region where the bulk of the 
7-ray emission is produced. 

In this work we compare Fermi /LAT 7-ray photon fluxes integrated over three different 
periods of time with Planck (sub-)mm observations of a sample of blazars. The data 
are presented in Section 2, and in Section 3 we perform a correlation analysis between 
intrinsic luminosities in both energy regimes. Section 4 presents an analysis of mm/sub-mm 
spectral shapes and 7-ray brightness. Our results are discussed and summarized in Section 
5. Throughout this manuscript, we use a ACDM cosmology with values within la of the 



WMAP results ( |Komatsu et al.||2009[ ); in particular, H =71 km s _1 Mpc^ 1 , Q m = 0.27, and 
tt A = 0.73. 



2. The data 



We build our analysis on the sample of 105 blazars presented in Giommi et al. 



(2011), who considered three samples of blazars with flux limits in the soft X-rays 



150 keV 



> 10 
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erg cm 



(count-rates 0.1-2.4 keV > 0.3 counts s _1 ), hard X-rays (S i 5 _ 
s _1 ), and 7-ray r] bands with additional 5 GHz flux density limits^] The advantage of 
using this sample is the availability of 7-ray flux measurements (quasi-) simultaneous to 
Planck observations. ESA's space mission Planck has been surveying the sky at nine 



frequencies since August 2009 (Planck Collaboration et al. 2011b). Its payload includes two 



instruments: the Low Frequency Instrument (LFI) operating at 30, 44 and 70 GHz while 
the High Frequency Instrument (HFI) observes at 100, 143, 217, 353, 545 and 857 GHz. 



1 selected from the Fermi/LAT Bright source list (Abdo et al. 



2009) with test statistics 



(TS) > 100 and Galactic latitudes larger than 10 degrees. 

2 soft X-rays, S 5GHz > 200 mJy; hard X-rays, $ 5G Hz > 100 mJy; 7-rays, S 5GHz > 1 Jy 
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The (sub-)mm flux densities employed in this work are listed in Table 6 of Giommi et al 



(2011) 



The Fermi /LAT 7-ray photon fluxes used in this study are based on data in the 
100 MeV to 300 GeV energy range and the data processing procedure has been fully 



described in section 3 of Giommi et al. (2011). For those sources where the Test Statistics 



(TS) in the whole energy band (0.1 - 300 GeV) delivers a value smaller than 25, an upper 
limit is given instead. Our gamma -ray data consist of three datasets which we describe 
below: 



• Simultaneous (< 5 7 > S im)'- The 7-ray photon flux was integrated over the period of 
time within which the source was observed at all Planck frequencies. A source in our 
sample typically sweeps over the Planck focal plane in about 2 weeks. Therefore, 
7-ray photon fluxes integrated during the Planck observation can be considered as 
simultaneous within one week. 

• Quasi-simultaneous (< S* 7 > qU a)'- The 7-ray photon flux densities were obtained by 
integrating Fermi data over a period covering two months and centered on the Planck 
observing period, i.e. about one month before and one month after the source was 
observed by Planck . We consider these 7-ray data as quasi- simultaneous to the 
Planck observations on monthly timescales. 

• Average (< S 1 > a ve)'- The 7-ray photon flux was integrated over a long-term period 
(27 months, from August 2008 to November 2010). We refer to these photon fluxes as 
average 7-ray fluxes. 

The number of 7-ray detections (Ndet) in each of the above Fermi datasets is as follows: 

< N det > sim = 26, < N det > qua = 56 and < N det > ave = 81. 
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Fig. 1. — Planck luminosities in (sub-) mm bands versus 27 months average Fermi lumi- 
nosities. Upper limits in the (sub-)mm and high-energy bands are indicated by left and 
downward arrows, respectively. Red circles represent FSRQs, green triangles BLLacs, and 
blue squares sources with uncertain types. The partial correlation parameters for each panel 
are given in Table [2} 

3. The (sub-)mm and 7-ray emission relationship 

In order to investigate the physical relation between the 7-ray and mm/sub-mm 
emission in blazars, we look for a statistical association between the wavelength regimes 
on the luminosity-luminosity plane. The reason for performing the correlation analyses 
using intrinsic luminosities rather than observed fluxes is that the flux-flux correlations 
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Fig. 2. — Planck luminosities in (sub-)mm bands versus 7-ray luminosities using Fermi data 
averaged over three different periods of time: 2 weeks (sim), 2 months (qua) and 27 months 
(ave). The partial correlation parameters for each panel are given in Table [3} 



will be distorted (and might even vanish) if the coefficient x in the relation L 7 oc L x v is 



different from unity (Feigelson & Berg 1983). On the other hand, if there is no intrinsic 



luminosity-luminosity correlation then no correlation will appear in the flux-flux plane as 



shown by Feigelson & Berg (1983). The mm/sub-mm luminosities have been computed 
according to the following expression, 



J (sub—)mm 



A 72 _1 



1 + * 



erg s 
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where is the luminosity distance, v corresponds to each of the nine Planck frequencies, 



and S u is the single-epoch flux density taken from Giommi et al. (2011). 



The 7-ray luminosity has been calculated according to 



(2) 



'1 + z) 1 -^ 

where ct 7 = V — 1, V being the 27 months average photon spectral index and 5*^,(0.1,300) 
is the photon energy flux integrated between 0.1 and 300 GeV calculated from the photon 
flux (Fry) by using the following expressions, 



£7(0.1,300) = 1.6 x 10 



1 — a, 



[3000 



1— a. 



1] a 7 7^ 1 



S 7 (0.1, 300) = 1.28 x lQ- d F 1 a 7 = 1 



(3) 
(4) 



derived from the relations presented in Ghisellini et al. (2009). 



3.1. Statistical methods 

It is conceivable that the dependence on redshift might induce artificial correlations 
between the luminosities. For that reason we must apply partial correlation methods to 
evaluate intrinsic L 7 — L u correlations. Because of the presence of upper limits in some of 
the 7-ray flux and (sub-)mm flux density measurements in our sample, survival analysis 
techniques are needed to properly quantify the correlation coefficient between these two 
emission bands. Therefore we have applied the Kendall's tau partial correlation for censored 



data (Akritas <__ Siebert 1996) to estimate whether there is an intrinsic correlation between 
luminosities once the influence of the upper limits has been taken into account and the 
common dependence on redshift has been removed. 



The FORTRAN code cens_taij^] implements the methods presented in Akritas 



3 http:/ /www2. astro, psu.edu/statcodes/cens_tau. f 



_9_ 



Siebert (1996) providing a measure of the correlation between the two luminosities while 
excluding the effects of the redshift (tllz)- The probability that non-correlation exists 
between luminosities (P T ) can be gleaned from the output of the code as well. We consider 
a partial correlation significant if the probability of non-correlation is less than 5% (P T < 
0.05). 

The slope of the correlations between (sub-)mm and 7-ray bands has been computed by 
the Akritas-Theil-Sen nonparametric regression method implemented in the routine cenken 
of the R package NADA. This method is best suited to our purposes owing the presence of 
upper limits in both variables. 



3.2. Significance test using the method of surrogate data 

To ensure that an intrinsic luminosity-luminosity correlation has not been drawn 
by chance we employed surrogate data (uncorrelated (sub-)mm and 7-rays luminosity 
pairs) to quantify whether the correlation coefficient obtained by applying the censored 
partial correlation test can or cannot be reproduced by mere coincidence. Based on the 



surrogate method described in detail in Ackermann et al. (2011a), we constructed surrogate 



luminosity data sets by following the next procedure: 

1. We start with a sample of m sources with a firm measurement of their redshift. The 
(sub-)mm and 7-ray luminosities are thus estimated by using expressions (1) and 
(2), respectively. Then, we run the censored partial correlation test to quantify the 
observed correlation strength (to). 

2. The sample is split in n redshift bins, each bin should contain at least 10 sources. This 



criterion follows the arguments described in section 4 of Ackermann et al. (2011a). 



3. In each bin we randomly shuffle the redshift and both the (sub-)mm and 7— ray 
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luminosities. The censored status of the luminosities at each band - either detected 
or upper limit - sticks to the original luminosity value. 

4. Using the surrogates from each bin, we build an uncorrelated data set of m pairs 
of luminosities with randomly assigned redshifts. Then we run the censored partial 
correlation test to estimate the correlation coefficient r» and the probability that the 
surrogate data set is non-correlated (Pj). 

5. We repeat the steps 3 and 4 a large number of times (N trials) and build the distribution 
of correlation coefficients for these surrogate uncorrelated luminosity pairs. 

We use the distribution of random correlation coefficients built on N trials = 1000 to 
compute the probability of obtaining a significant correlation (Pj < 0.05) with a coefficient 
larger or equal (in terms of absolute value) than the observed one. Then, the significance 
of the observed censored partial correlation (P surrogate) is estimated by computing the ratio 
N hits /N trials, where Nuts is he number of times that we registered a significant correlation 
with |Tj| > Tq. 

An observed censored partial correlation is considered real only if the probability of 
getting the same (or a larger) correlation coefficient is less than 5% {P surrogate < 0.05). This 
criterion will ensure us that observing a significant (sub-)mm ad 7-ray luminosity-luminosity 
partial correlation with a coefficient tq is not likely to occur by chance. 



3.3. The L 7 - L( SU 5_) mm correlation 



We investigate the luminosity- luminosity correlation between the (sub-)mm emission 
and 7-ray fluxes in the sample of blazars presented in Giommi et al. (2011). We only 
consider sources with firm measurements of their redshift and 7-ray photon flux - either 
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detection or upper limit- leaving us with a total of 98 sources. The luminosities at the 
(sub-)mm and 7-ray regimes were computed according to expression 1 and 2, respectively. 
In equation (2), Sl^O.l, 300) is the photon energy flux integrated over 27 months (< S ave >). 
For those cases where non-detection of the sources was possible, we estimate an upper limit 
of the luminosity assuming the average photon spectral index for our sample, < T >= 1.38. 
All luminosities are listed in Table [U 

Figure 1 shows the luminosity-luminosity relation between all (sub-)mm frequency 
bands and average values of 7-ray emission obtained by integrating Fermi data over 27 
months. In order to examine the possible dependence of the correlation on blazar type, the 
various blazar types are symbol coded as shown in the legend of Figure 1. FSRQs populate 
the upper right part of the luminosity- luminosity plane, whereas BLLacs and sources with 
uncertain type continue the trend toward lower luminosities. It should be noticed that some 
of the sources classified as "uncertain type" are actually known radio galaxies, however, for 



sake of consistency we keep the nomenclature as stated in Tables 1 to 3 of Giommi et al 



(2011) 



The correlation parameters and their significances for each of the panels of Figure 1 
are summarized in Table El The statistical methods described in the above subsections 
reveal that mm/sub-mm luminosities are positively correlated with 7— ray luminosities over 
five orders of magnitude when all source types are considered. However, differences arise 
when we compute partial correlations for FSRQs and BLLacs, separately. The correlations 
between (sub-)mm and 7-ray emission bands are significant and real when FSRQs are 
considered only. However, the surrogate test method does not provide evidence that 
correlations for BLLacs are significant. 

The slope of the relation between (sub-)mm and 7-ray luminosities has been computed 
with the Akritas-Theil-Sen nonparametric regression method and the fitted values are 
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also listed in Table 1. For all 98 sources considered in our study, the (sub-)mm - 7- ray 
luminosity relation can be well approximated as L 7 ~ L*, where < xau >— 1.27 ±0.02. The 
slope fitted for BLLacs (< xblluc >= 1.13 ± 0.03) is somewhat shallower than the slope 
computed for FSRQ (< xfsrq >= 1-40 ± 0.03). While the different photon spectral index 



between FSRQs and BLLacs (Ghisellini et al. 2009) may play some role in the difference 



between fitted slopes, the lack of correlation for BLLacs might indicate that different 
radiation mechanisms are behind the 7-ray emission. As we discuss in section 4, clean and 
well characterized samples of BLLacs are needed to get further insight into the physical 
conditions that make the difference between FSRQs and BLLacs. 

The low detection rate at the highest Planck frequencies can be a combination of the 



high flux density detection limit at these frequencies - see Table 3 in Planck Collaboration 



et al. (2011c), and the fact that blazars in general become brighter at sub mm frequencies 



only during flaring states (e.g Marscher & Gear 1985 Raiteri et al. 2011). Therefore, due 



to the relatively low number of sources detected at 545 and 857 GHz the correlation and 
best-fitted line parameters for these frequency bands shown in Table [2] and [3] should be 
taken with caution. 



3.4. Dependence of the L 7 - L( su b_) mm correlation on simultaneity 

We find that in our sample of blazars the 7-ray emission averaged over large periods 
of time (i.e. 27 months) is significantly correlated to the (sub-)mm radiation. This result 
suggests a coupling of the emitting regions, however, to get a physical insight about this 
relation we next need to investigate whether the correlation between (sub-)mm and 7-ray 
emission bands shows a trend for strengthening with simultaneity. We emphasize again 
that our Fermi data set contains information about the 7-ray emission averaged over three 
different periods of time, allowing us to assess the intrinsic correlation between (sub-)mm 
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radiation and simultaneous, quasi-simultaneous and averaged 7-ray emission. 

To investigate the dependence of the L 1 - L( su b_) mm correlation on simultaneity of the 
data, the levels of censoring at 7-rays must be controlled and therefore we only consider 
sources that were detected by Fermi at each of the three averaging periods. Hereafter we 
refer to this subsample as the Fermi- detected sample which comprises 24 sources. Most of 
them are classified as FSRQ type and based on their simultaneous SEDs all of them have 



been classified as low-synchrotron peaked (LSP) in Giommi et al. (2011). 



Figure 2 shows the luminosity-luminosity relation for all (sub-)mm frequency bands 
and 7-ray emission integrated over three different periods of time. The simultaneity of the 
7-ray emission to the Planck measurements is symbol coded as shown in the legend at the 
right center of the plot. Table [3] reports the partial correlation parameters and associated 
probabilities for the three 7-ray data sets. 

As can be seen from Table |3j it is noticeable that the strongest correlations arise when 
Fermi measurements have been integrated over long periods of time (2 and 27 months), 
whereas the weakest and less significant correlations are always obtained when simultaneous 
observations are involved. The remarkably stronger correlation with the 7-ray photon fluxes 
integrated over 27 months is likely an artificial effect due to averaging the 7-ray photon 
fluxes over a very long period of time. This reduces the dynamical range of 7-ray photon 



flux observed (Muecke et al. 1997), which in turn reduces the scatter of the L^-L v relation 



and improves the correlation strength. Given the evidence that the strongest 7-ray flares 



tend to have a short duration time, with a typical timescale of a month (see Abdo et al. 



2010b), one might argue that by averaging 7-ray data over two months we may be averaging 
flares, which would reduce the dynamical range of quasi-simultaneous 7-ray luminosities, 
and therefore stronger correlations may be induced. However, it is not likely that all sources 
considered in our study were flaring during the weeks that Planck observed them. This 



-14- 



leads us to believe that the the most significant correlations in our study arise when using 
the 7-ray photon fluxes averaged over 2 months. 



4. Shape of the synchrotron spectra and 7-ray brightness 

In this section we explore emerging trends between (sub-) mm spectral shape and 
7-ray brightness with the aim to find out whether the (sub-) mm spectrum shape can be 
considered as a new piece to the puzzle of what makes a blazar shine at 7-rays. The radio to 



sub-mm spectra can be well approximated by two power laws (e.g. Gear et al. 1994; Planck 



Collaboration et al. 2011a), however the point where both power laws connect has always 
been assumed rather arbitrarily based on the spectral coverage of the data sets involved. 
This in turn might introduce some biases to the spectral indices fitted to parametrize the 
spectral shapes. 

Hence, in this work we use a different approach by fitting the Planck observations 
from 30 to 857 GHz with a broken power law model where the fiducial point where the 
intersection of the power laws occurs is considered as a free parameter along the fit. The 
expression of the broken power law used to model the (sub-)mm spectral range is given by, 



V amm for V < U b reak 

S(u) oc { (5) 

V^Zk asubmm) V a ™ b ™ for V > V break , 

where a mm and a su b mm are the spectral indices for the mm and sub-mm part of the 
spectrum, respectively, and the two power laws are connected at the break frequency z/ breafc . 

To perform the analysis of spectral shapes and 7-ray brightness, we select only those 



sources from Giommi et al. (2011) that have been detected in at least five Planck frequency 



bands, this with the aim to produce a reliable modeling of the spectra. The 47 sources 
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selected based on the criterion mentioned above, along with the best-fit model parameters, 
are listed in Table HI 

The new (sub-)mm observations provided by the Planck satellite and their modeling 
with a broken power law function allow us to typecast the (sub-)mm spectral shapes in 
our sample. Figure [3] shows all four spectral shapes that could possibly be found in our 
sample. Next we describe their characteristics and possible interpretation in terms of 
multicomponent synchrotron spectra: 

(A) Steep a mm - Steep a su b mm : (sub-)mm emission decreases monotonically with 
frequency. This straight spectral shape from mm to sub-mm regimes is dominated 
by the underlying emission from a component with a synchrotron turnover frequency 
located at cm wavelengths. Quiescent activity levels (i.e. no new young synchrotron 
components) is expected for sources featuring this spectral type. 

(B) Steep a mm - Flat/rising a su b m m : mm spectrum is dominated by the emission from 
an aged component that becomes self absorbed at the cm regime, as in (A). Despite 
that, a flat/rising a su b mm , spectrum betrays the presence of a new component in an 
early development stage, located well within the radio core region. Then a new flare 
connected to the ejection of a new VLBI component is likely to occur before the 
sub-mm spectrum steepens. It is possible that a rising a su b mm could also be related to 
the presence of a dust component, however, a more detailed dissection of the spectra is 
needed to assess its presence. 

(C) Flat a mm - Steep a su b mm : An overall flatness of the mm spectrum can be explained 
by a succession of several components in the jet as stated in the so-called "cosmic 



conspiracy" (Cotton et al. 1980). Nevertheless, the steepening of the sub mm spectra 
is consistent with synchrotron emission from a single component which becomes self 
absorbed somewhere in the middle of the mm wavelength regime. 
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(D) Rising a mm - Steep a su b m m : In the absence of adjacent strong spectral components 
at lower or higher frequencies, a new component recently ejected from the radio core, 
with a self absorbed spectrum peaking in the mm regime, dominates the overall shape 
of the radio spectrum. 

The (sub-)mm spectra of the selected blazars were classified into the above spectral 
categories based on the relative difference between spectral indices Aa = " s "'" nm ~" mm . The 

Otrrvm 

spectral classification for the 47 sources selected is reported in Table [4] and is defined by the 
following conditions: 

• A : a mm ,a su b mm < and |Aa| < 0.5; relative difference between indices is less than 
50%, 

• B: a mm < and a submm > 0, 

• C : a mm , a su bmm < and |Aa| > 0.5; relative difference between indices is greater 
than 50%, 

• D: a mm > and a submm < 0, 

The classification scheme described above is sketched in Figure |3j which attempts to 
compare the (sub-)mm spectral shapes of blazars at different stages of the synchrotron 
components evolution. Based on successive approximations, we further group the four 
categories into two main spectral classes: (i) A-B includes sources previously classified as 
A or B , these represent the quiescent and the very initial stages of an ejected component, 
respectively, (ii) The spectral type C-D is conformed by sources whose spectral shapes are 
categorized as type C or D, these specific spectral features indicate the presence of a freshly 
ejected synchrotron component where the spectral turnover is somewhere in the middle of 
the mm wave regime. 
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Based on the latter spectral type classification, we may classify our sample into 23 
sources of category A-B and 24 sources falling into the category of C-D. Figures [4] and [5] 
show the (sub-)mm spectra for both spectral types and overplotted are the best-fits to the 
(sub-)mm spectra using the broken power law model. 

Due to the fact that L 7 — L v are best correlated when using quasi-simultaneous 
7-ray fluxes (< 5* 7 > qua ), we proceeded to find out if there is a pattern in the spectral 
shapes that relates to the levels of quasi-simultaneous 7-ray emission. Figure [6] shows the 
quasi-simultaneous 7-ray flux distributions as boxplots. The levels of 7-ray emission for 
sources with (sub-)mm spectra belonging to the spectral classes A and B (spectral type 
A-B) are displayed in the left boxplot of Figure [6] while sources with spectra of the type 
C and D (spectral type C-D) are shown in the boxplot to the right. Upper limits and 
detections are symbol coded as denoted in the legend of Figure [6j 

At first glance, it seems that sources of the spectral type C-D are brighter at 7-rays 
than sources of the spectral type A-B. This can be gleaned from the size of the boxplots 
and the locations of the medians represented by the thick line at the middle of the boxes 
displayed in Figure |6j Differences in the distribution of quasi-simultaneous 7-ray photon 
fluxes for sources with different (sub-)mm spectral types were statistically investigated. 
Because of the presence of upper limits to the 7-ray photon fluxes we apply the univariate 



two sample tests implemented in the ASURV package (Isobe et al. 1986; Lavalley et al 



1992) to the 7-rays photon flux distributions displayed in Figure |6| We further consider that 
two populations are significantly different if the probability that the differences between 
populations arise by chance is P < 0.05. 

All five two-sample tests implemented in ASURV (Gehan's generalized Wilcoxon test 
- permutation variance, Gehan's generalized Wilcoxon test - hypergeometric variance, 
logrank test, Peto-Peto generalized Wilcoxon test and Peto-Prentice generalized Wilcoxon 
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> 

Frequency (GHz) 



Fig. 3. — An illustration of all synchrotron spectral shapes we could possibly find in our 
sample, see text for detailed description. The overall spectrum has been approximated 
by a broken power law with spectral indices a mm (red) and a su b- mrn (blue), respectively. 
The break frequency (inbreak) is indicated by the grey thick vertical line. Aged synchrotron 
components are represented by dashed lines 
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test) converge to the same result: the levels of 7-ray emission of spectral types A-B and 
C-D are significantly different and the probability that such difference arises by chance is 
less than 4%. The statistical tests show that the shapes and medians of the distributions 
of 7-rays for spectral types A-B and C-D are significantly different, where the levels of 
7-ray emission for sources belonging to the type C-D are significantly higher than for those 
classified as A-B. 

We note that other classification schemes for blazars have already been proposed based 



on the variability of the radio continuum spectra from cm to mm wavelengths (Angelakis 



et al. 2011). However, the classification scheme proposed in this work - based on successive 
approximations of the synchrotron components evolution - for the fist time permits to use 
the (sub-)mm spectral shapes as flags for high levels of 7-ray emission. 



5. Discussion and summary 

While various studies during the Fermi-era have found evidence for a correlation 



Kovalev et al. 


2009 


Ackermann et al. 


2011a 


Ghirlanda 



et al. 2011 Nieppola et al. 2011), our work addresses for the first time this connection at 
(sub-)mm wavelengths. We have found a significant correlation between (sub-)mm and 7-ray 
luminosities over five orders of magnitude. Since we have made use of partial-correlation, 
surrogate data methods and survival analysis techniques we are confident that the 
correlation L 1 oc L v is not an artifact of the detection limits, chance or due to the common 



dependence on redshift. This result is consistent with previous studies (e.g. Ackermann 



et al. 2011a; Ghirlanda et al. 2011 Nieppola et al. 2011 Arshakian et al. 2012), albeit 



using high frequency radio observations. Giommi et al. (2011) do not report evidence of 



a significant correlation between 143 GHz and 7-ray fluxes in some of the sub-samples 
considered. This apparently contradiction with our findings is not likely to be influenced 
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Fig. 4. — The (sub-)mm spectra for 23 sources classified within the spectral type A-B. See 
Figure [3] for a visual description of the types. 



-21- 




10 100 1000 10 100 1000 10 100 1000 10 100 1000 



v (GHz) 



Fig. 5. — The (sub-)mm spectra for 24 sources classified within the spectral type C-D. See 
Figure [3] for a visual description of each spectral type. 
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Fig. 6. — The distribution of of 7-ray emission levels for sources with spectral types A- 
B (left) and C-D (right) are shown in boxplots, where the thick line at the middle of the 
box represents the median. The boxes comprise 75% of the distributions. The levels of 
7-ray emission for sources belonging to the type C-D are significantly higher than for those 
classified as A-B 
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by sources without redshift not considered in this study. Instead, such difference may be 
attributable to the fact that correlations in the flux-flux plane can vanish if the slope of the 
luminosity-luminosity correlation is different from unity (see discussion in Section 3) as it 
turns out to be the case for the L 7 oc L^ sub _) mm relation (see Table |2]). 

The positive correlation between luminosities is still robust and becomes even stronger 
when only FSRQs are considered, however, the correlation vanishes for BLLacs as a group. 
In addition, the slope of the fitted L 7 cx L x v relation is shallower for BLLacs. The different 
gamma-radio behavior between blazar types has been found in previous studies using 



37 GHz data (Leon-Tavares et al. 2011 Nieppola et al. 2011). On the other hand, previous 



studies at low frequencies (Ghirlanda et al. 2011 Ackermann et al. 2011a) find that the 



radio-gamma correlation is significant for both FSRQs and BLLacs, separately. This 
apparent discrepancy might be related to the fact that submm spectral indices of FSRQs 



and BLLacs have been found to be significantly different (Gear et al. 1994). The latter 



result was interpreted as an intrinsic difference in the underlying jets of BLLacs and FSRQ. 
Unfortunately due to the low number of BLLacs in the sample analyzed in section 4 - see 
Table |4j we could not explore this difference with the new (sub-)mm observation. 



Nevertheless, before speculating into real physical differences between FSRQs and 
BLLacs, we would like to rise a word of caution regarding the correlation results solely for 



the BLLac population. Giommi et al. (2012) have drawn attention to the fact that some 



sources classified as BLLacs might have been misclassified with FSRQ. This issue can be 
a direct consequence of the rather artificial criterion that draws the line between FSRQs 
and BLLacs - rest frame equivalent width (EW) less than 5A - introduced twenty years 



ago (Stickel et al. 1991). It is plausible that a particular source was originally classified 



as BLLac at the time when its non-thermal optical continuum was prominent and high 
enough to swamp any broad emission line in the optical spectral range, as recently shown 
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in Padovani et al. (2012). However, after some time the same source could have been 



classified as FSRQ when the optical continuum emission weakened, allowing the exhibition 
of broad emission lines. To overcome this misclassification problem we should adopt a 
new blazar taxonomical criterion to classify sources with intrinsically prominent and weak 



broad-line regions as recently proposed in Ghisellini et al. (2011) and Giommi et al. (2012). 



Therefore a compilation of a large sample of true BLLacs whose parent population includes 



low-excitation radio galaxies (LERG) displaying FR-I morphologies, as proposed in Giommi 



et al. (2012), is highly desirable to draw robust conclusions on the (sub-)mm and 7-ray 



connection in BLLacs. 

We now turn to the dependence of the L 7 oc L v relation on simultaneity. Although 
the strongest correlation parameters are achieved when using Fermi 7-ray photon fluxes 
averaged over 27 months, we believe that such tight correlations are a consequence of 
averaging the 7-ray photon fluxes over a very long period of time. This reduces the scatter 
of the L^-Ly relation while the detection of faint 7-ray sources increases the 7-ray luminosity 
dynamical range, hence the correlation strengths are artificially enhanced. Furthermore, the 
L 1 oc L u relation becomes weak if 7-ray and (sub-)mm data are obtained simultaneously. 
On the other hand, including quasi-simultaneous observations within monthly timescales 
produces a more significant correlation. This might suggest an intrinsic delay between 
the emission in both wavebands, however, with the data currently available we cannot 
disentangle whether the (sub-)mm emission leads the production of 7-rays or vice versa. 
There is, however, recent evidence of delays on monthly scales between mm and 7-ray 
emission bands, the onset of a mm-flare occurring about 2 months before the maximum of 



the 7-ray emission is reached ( Leon-Tavares et al. 2011). 



Perhaps the most remarkable result we find is that brighter 7-ray sources tend to show 
a specific shape in their radio spectra, i.e. flat/rising mm and falling sub-mm spectral 



- 25 - 



indices (spectral type C-D), see Figure [3] for the classification scheme proposed in this 
work. This spectral shape is consistent with a single synchrotron component becoming self 
absorbed in the middle of the millimeter wavelength regime that must be responsible for 
the last spectral turnover. Despite the fact that we have not attempted to properly dissect 
the (sub-)mm spectra into single synchrotron components, turnover frequencies for sources 
of the spectral type C-D can be roughly constrained to u m > 30 GHz. 

Our finding that 7-ray emission is enhanced in sources of the spectral type C-D - 
featuring an approximate spectral turnover at the middle of the mm wavelength regime 
- suggests an important role of the plasma jet disturbances, responsible for shaping the 
(sub-)mm spectra, in producing the 7-ray emission. Since B oc ( Lobanovj|l998 ) sources of 



the spectral type C-D [y m > 30 GHz) are expected to harbor strong magnetic fields in their 
cores leading to polarized emission. This argument finds support in recent observational 



results by Linford et al. (2012), where it was found that strong core polarization is related 



to intense 7-ray emission. 

The existence of large magnetic fields may reflect in either large electron densities in 
the emitting region (in this case spectra must be dominated by emission from the core) 
or the presence of growing shocks emerging from the radio core region. To disentangle 
the nature of the emission shown in the snapshot (sub-)mm spectra of Figures [i] and [5j 
multifrequency VLBI monitoring is highly desirable to properly dissect the spectra into 



individual components (e.g. Lobanov h Zensus 1999 Savolainen et al. 2008) and to identify 



the contribution of the core and the emerging shocks to the overall radio spectrum. In 
addition, upcoming Planck data releases will allow us to build multi-epoch (sub-)mm 
spectra. Then, by studying the evolution of the spectral shapes, we can get an insight 
into the changes of the physical properties of the jet (e.g., density, magnetic field, Doppler 
factor), and those must serve as input to theoretical models aiming to reproduce the levels 
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of 7-ray emission observed. 

In summary, we have used a joint sample of soft X-ray, hard X-ray and 7-ray selected 
blazars to explore for the first time the relationship between (sub-)mm emission from 30 to 
857 GHz and 7-ray photon fluxes integrated over three different periods of time: when a 
source was in the Planck field of view, one month before and after the source was observed 
by Planck, and over 27 months. Our findings can be listed as follows: 

• We have found a significant correlation between 7-ray and (sub-)mm luminosities 
which holds over five orders of magnitude. Moreover, using a sample of bright 7-ray 
sources detected by Fermi on each of the averaging periods mentioned above, we 
find - based on our statistical analysis and bearing in mind that 7-ray photon fluxes 
averaged over long periods of time (i.e. 27 months) enhance artificially the correlation 
strength - that quasi- simultaneous observations (within 2 months) showed the most 
significant statistical correlations. 

• Sources with high levels of 7-ray emission show a characteristic signature in their 
radio spectra: flat/rising mm and steep sub-mm spectral indices (spectral type C-D in 
Figure [3]). This spectral shape can be associated with a single synchrotron component 
which becomes self absorbed in the middle of the mm wavelength regime. Such high 
spectral turnover frequencies reveal the presence of emerging disturbances in the jet 
which are likely to be responsible for the high levels of 7-ray emission. 

We thank the anonymous referee for careful reading and helpful comments. The 
Metsahovi team acknowledges the support from the Academy of Finland (project numbers 
212656, 210338, 121148, and others). We acknowledge the use of data and software facilities 
from the ASI Science Data Center (ASDC), managed by the Italian Space Agency (ASI). 
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Table 2. Results of the correlation analysis between the mm/sub-mm and average 7-rays 

luminosities 



< Ly >ave 



n/v-lp lanck /ul FFrmi T llz P T Piurrojote Z 

(1) (2) (3) (4) (5) 



L 30 GHz 



ALL 


98/31/23 


0.38 


< 10~ 6 


< lO" 4 


1.26 


FSRQ 


51/10/10 


0.40 


5.X10 -5 


< 10~ 4 


1.45 


BLLac 


27/14/3 


0.34 


9X10 -5 


0.1 


1.09 


'44 GHz 












ALL 


98/46/23 


0.29 


2xl0~ 7 


< 10- 4 


1.30 


FSRQ 


51/18/10 


0.33 


5xl0~ 4 


lxl0~ 3 


1.46 


BLLac 


27/18/3 


0.23 


7x 10~ 4 


0.1 


1.18 


'70 GHz 












ALL 


98/43/23 


0.30 


6xl0~ 8 


< 10- 4 


1.29 


FSRQ 


51/16/10 


0.35 


5X10 -5 


< 10~ 4 


1.46 


BLLac 


27/18/3 


0.24 


4x 10 -3 


0.1 


1.17 


■inn GHz 












ALL 


98/33/23 


0.31 


2xl0~ 6 


< 10- 4 


1.21 


FSRQ 


51/13/10 


0.33 


7xl0~ 4 


2xl0~ 3 


1.30 


BLLac 


27/13/3 


0.25 


7x 10 -3 


0.2 


1.06 


'143 GHz 












ALL 


98/31/23 


0.34 


lxl0~ 6 


< lO" 4 


1.19 


FSRQ 


51/11/10 


0.37 


2xl0~ 4 


< lO" 4 


1.30 


BLLac 


27/14/3 


0.28 


BX10 -3 


0.1 


1.03 


'217 GHz 












ALL 


98/33/23 


0.32 


4xl0~ 7 


< 10- 4 


1.21 


FSRQ 


51/13/10 


0.34 


3xl0~ 4 


1X10 -3 


1.31 


BLLac 


27/12/3 


0.28 


3X10 -3 


0.1 


1.04 


'353 GHz 












ALL 


98/51/23 


0.21 


lxl0~ 5 


< 10- 4 


1.31 


FSRQ 


51/25/10 


0.24 


2xl0~ 3 


2xl0~ 3 


1.41 


BLLac 


27/15/3 


0.25 


lxl0~ 3 


0.1 


1.13 


'545 GHz 












ALL 


98/70/23 


0.13 


8xl0~ 4 


< 10" 4 


1.36 


FSRQ 


51/33/10 


0.17 


0.01 


6xl0~ 3 


1.5 


BLLac 


27/22/3 


0.11 


0.05 


0.2 


1.23 


'857 GHz 












ALL 


98/83/23 


0.10 


3xl0~ 3 


< 10~ 4 


1.32 


FSRQ 


51/42/10 


0.11 


5xl0~ 2 


lxl0~ 3 


1.36 


BLLac 


27/23/3 


0.08 


8xl0~ 2 


0.1 


1.24 



Note. — Column (1): (number of sources / number of Planck upper limits / 
number of Fermi upper limits ). Column (2) th z : partial correlation coefficient. 
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Column (3) P T : probability that non-correlation exists between luminosities after 
the common effect of redshift has been removed . Colum (4) P S urrogate'- probability 
that the correlation arises by chance (estimated by surrogate methods). Column (5): 
slope fitted to the luminosity-luminosity relation. 



-37- 



Table 3. Results of the correlation analysis between the mm/sub-mm and 7-rays emission 

integrated over different periods of time 







< > s i m 




< 


> qua 






< 




m,z 


P-r P surrogate 


T ll,z 


Pr 


P surro^ 


• ate 


m,z 


P T 


P surrogate 


L30 


0.45 


0.03 3xl0~ 3 


0.48 


0.02 


2x10" 


-3 


0.50 


0.01 


< lO" 4 


L44 


0.34 


0.06 


0.35 


0.04 


6x10" 


-3 


0.39 


0.03 


lXlO -3 


L70 


0.46 


0.02 3X10 -3 


0.46 


0.01 


< 10- 


■4 


0.47 


0.01 


< 10- 4 


Lino 


0.36 


0.06 


0.39 


0.03 


lx 10" 


-2 


0.39 


0.03 


lXlO" 2 


L143 


0.38 


0.06 


0.41 


0.04 


7x10" 


-3 


0.46 


0.02 


< 10- 4 


L217 


0.34 


0.08 


0.36 


0.06 






0.39 


0.04 


5X10 -3 


L353 


0.32 


0.08 


0.34 


0.06 






0.38 


0.04 


2xl0~ 3 


L545 


0.32 


0.04 < 10~ 4 


0.33 


0.03 


< 10" 


■4 


0.35 


0.03 


< 10~ 4 


L8S7 


0.29 


0.04 < 10 -4 


0.29 


0.03 


1x10" 


-3 


0.26 


0.05 


< lO" 4 


Note. 




partial correlation coe 


fficient. P T : 


probability that 


non- 


-correlation 


exists between luminosi- 



ties after the common effect of rcdshift has been removed . P S urrogate- probability that the correlation arises 
by chance (estimated by surrogate methods). 
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Table 4. (Sub-)mm spectral fit parameters 



,12000 Name ctmm, a subm,m ^break Spectral class 

[GHz] 



0010+1058 


-0.61 


0. 


58 


343.13 


B 


0120-2701 


-0.09 


-0 


.03 


100.00 


B 


0136+4751 


-0.24 


-0 


.52 


120.57 


C 


0210-5101 


-0.06 


-0 


.86 


83.40 


C 


0217+0144 


-0.08 


-0 


.40 


142.82 


c 


0217+7349 


-0.66 


-0 


.74 


70.06 


A 


0221+3556 


-0.52 


0. 


24 


159.85 


B 


0237+2848 


-0.35 


-0 


.88 


128.39 


C 


0238+1636 


0.21 


-0 


.60 


111.08 


D 


0319+4130 


-0.50 


-0 


.76 


138.77 


C 


0336+3218 


-0.18 


-0 


.56 


91.92 


C 


0423-0120 


-0.25 


-0 


.61 


197.82 


c 


0428-3756 


-0.43 


-0 


.63 


77.46 


A 


0457-2324 


-0.37 


-0 


40 


198.48 


B 


0522-3627 


-0.09 


-0 


.42 


110.00 


C 


0530+1331 


0.47 


-0 


.75 


59.04 


D 


0538-4405 


-0.02 


-0 


.32 


130.00 


C 


0539-2839 


-0.53 


-0 


.65 


142.97 


A 


0841+7053 


-0.11 


-0 


.97 


96.84 


D 


0854+2006 


-0.05 


-0 


52 


81.45 


B 


0920+4441 


-0.39 


-0 


.09 


303.81 


B 


1058+0133 


-0.44 


-0 


.52 


217.81 


A 


1058-8003 


-0.14 


-1 


12 


95.28 


C 


1130-1449 


-0.44 


-0 


.33 


180.09 


B 


1229+0203 


-0.13 


-0 


.96 


126.62 


C 


1246-2547 


0.34 


-0 


.46 


83.23 


D 


1256-0547 


-0.16 


-0 


.56 


107.54 


C 


1310+3220 


-0.63 


-0 


61 


195.47 


A 


1504+1029 


-0.37 


-0 


.97 


87.73 


C 


1517-2422 


-0.24 


0. 


04 


180.06 


B 


1751+0939 


-0.46 


-0 


.57 


127.79 


A 


1800+7828 


-0.30 


-0 


.35 


183.88 


A 


1833-2103 


-0.42 


-0 


.61 


111.17 


A 


1911-2006 


-0.40 


0. 


39 


296.73 


B 


1923-2104 


0.11 


-0 


48 


78.37 


D 


1924-2914 


-0.42 


-0 


.60 


90.06 


A 


2147+0929 


-0.17 


0. 


00 


100.00 


B 


2148+0657 


-0.59 


-0 


.68 


205.52 


A 


2151-3027 


-0.03 


-0 


.57 


99.98 


C 


2202+4216 


-0.21 


-0 


.64 


113.32 


C 


2203+3145 


-0.08 


-0 


.76 


109.95 


C 


2207-5346 


-0.37 


-0 


.41 


124.33 


B 


2229-0832 


0.24 


-0 


.64 


96.46 


D 


2232+1143 


-0.56 


-0 


.92 


99.93 


C 


2253+1608 


0.47 


-0 


.26 


97.99 


D 
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Table 4 — Continued 



J2000 Name 




a subTnm 


"break 


Spectral class 








[GHz] 




2327+0940 


-0.70 


-0.82 


129.50 


A 


2333-2343 


-0.13 


-0.09 


336.43 


B 



Note. — The (sub-)mm spectra have been approximated with 
the broken power law model of equation (5). 



